****************************************************************************
* ``The Costs and Environmental Justice Concerns of NIMBY in Solid Waste Disposal''*
* Phuong Ho (SNF-Centre for Applied Research in Norwegian School of Economics) *
*******************************************************************************
*** This file makes Figure 2 (goodness of fit) ***************************
clear all 
set more off

global mydir "C:\XYZ" // set working directory
global data "$mydir\data"
global matlabin "$mydir\matlab\input"
global result "$mydir\TEX\"
cd "$mydir\stata"
set matsize 11000
set maxvar 32700

***********************************************************************
**** IMPORT STATUS QUO ESTIMATES *************************************
**********************************************************************
import excel "$mydir\matlab\output\nimby_statusquo.xlsx", sheet("Sheet1") firstrow clear
label var sj_q "market size"
label var pr "model implied market share"
label var Bq "model implied trash amount"
label var el "price elasticity"
label var CS "hauler's surplus"
sort TCJ
***** GOODNESS OF FIT ***********************************************
***** FIGURE 3 ******************************************************
// Subfigures A and B. Correlation in trash quantity
gen qInThousand=q/1000
gen predictedQInThousand=Bq/1000
label var qInThousand "data flow"
label var predictedQInThousand "predicted flow"
correlate qInThousand predictedQInThousand
twoway scatter predictedQInThousand qInThousand || line qInThousand qInThousand, sort ///
   graphregion(color(white)) legend(off) ///
   xtitle("observed flow (thousand tons)",placement("east")) ///
   ytitle("predicted flow (thousand tons)", placement("north"))
graph export "$result\figure2a.pdf", replace

correlate qInThousand predictedQInThousand if year==2010
twoway scatter predictedQInThousand qInThousand if year==2010 || line qInThousand qInThousand  if year==2010, sort ///
   graphregion(color(white)) legend(off) ///
   xtitle("observed flow (thousand tons)",placement("east")) ///
   ytitle("predicted flow (thousand tons)", placement("north"))
graph export "$result\figure2b.pdf", replace

gen distanceQWeighted=q*dridis
gen distanceQWeighted_predicted=Bq*dridis
gen priceQWeighted=q*realprice
gen priceQWeighted_predicted=Bq*realprice
sort year quarter county Code
by year quarter county: egen tdistanceQWeighted=total(distanceQWeighted)
by year quarter county: egen tq=total(q)
by year quarter county: egen tdistanceQWeighted_predicted=total(distanceQWeighted_predicted)
by year quarter county: egen tqf=total(Bq)
by year quarter county: egen tpriceQWeighted=total(priceQWeighted)
by year quarter county: egen tpriceQWeighted_predicted=total(priceQWeighted_predicted)
by year quarter county: gen dup=cond(_N==1,0,_n)

replace tdistanceQWeighted=tdistanceQWeighted/tq
replace tdistanceQWeighted_predicted=tdistanceQWeighted_predicted/tqf
replace tpriceQWeighted=tpriceQWeighted/tq
replace tpriceQWeighted_predicted=tpriceQWeighted_predicted/tq

// Subfigure C and D. Correlation in waste-weighted distance
correlate tdistanceQWeighted tdistanceQWeighted_predicted if dup<=1
twoway scatter  tdistanceQWeighted_predicted tdistanceQWeighted if dup<=1 || line tdistanceQWeighted tdistanceQWeighted if dup<=1, sort ///
   graphregion(color(white)) legend(off) ///
   xtitle("observed distance (mile)",placement("east")) ///
   ytitle("predicted distance (mile)", placement("north"))
graph export "$result\figure2c.pdf", replace

correlate tdistanceQWeighted tdistanceQWeighted_predicted if year==2010  & dup<=1
twoway scatter tdistanceQWeighted_predicted tdistanceQWeighted if year==2010  & dup<=1 || ///
	line tdistanceQWeighted tdistanceQWeighted if year==2010  & dup<=1, sort ///
   graphregion(color(white)) legend(off) ///
   xtitle("observed distance (mile)",placement("east")) ytitle("predicted distance (mile)", placement("north"))
graph export "$result\figure2d.pdf", replace

// Subfigure E and F. Correlation in waste-weighted price
correlate tpriceQWeighted tpriceQWeighted_predicted if dup<=1
twoway scatter  tpriceQWeighted_predicted tpriceQWeighted if dup<=1 || line tpriceQWeighted tpriceQWeighted if dup<=1, sort ///
   graphregion(color(white)) legend(off) ///
   xtitle("observed price ($/ton)",placement("east")) ///
   ytitle("predicted price($/ton)", placement("north"))
graph export "$result\figure2e.pdf", replace

correlate tpriceQWeighted tpriceQWeighted_predicted if year==2010  & dup<=1
twoway scatter tpriceQWeighted_predicted tpriceQWeighted if year==2010  & dup<=1 || ///
	line tpriceQWeighted tpriceQWeighted if year==2010  & dup<=1, sort ///
   graphregion(color(white))  legend(off) ///
   xtitle("observed price ($/ton)",placement("east")) ///
   ytitle("predicted price ($/ton)", placement("north"))
graph export "$result\figure2f.pdf", replace
